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Abstract 



One of the challenges in epidemiology is to account for the complex morphological structure of hosts such 
Oh ■ as plant roots, crop fields, farms, cells, animal habitats and social networks, when the transmission of 

Q I infection occurs between contiguous hosts. Morphological complexity brings an inherent heterogeneity 

in populations and affects the dynamics of pathogen spread in such systems. We have analysed the 
influence of realistically complex host morphology on the threshold for invasion and epidemic outbreak 
^ in an SIR (susceptible-infected-recovered) epidemiological model. We show that disorder expressed in 

the host morphology and anisotropy reduces the probability of epidemic outbreak and thus makes the 
system more resistant to epidemic outbreaks. We obtain general analytical estimates for minimally safe 
^ I bounds for an invasion threshold and then illustrate their validity by considering an example of host 

data for branching hosts (salamander retinal ganglion cells). Several spatial arrangements of hosts with 
different degrees of heterogeneity have been considered in order to analyse separately the role of shape 

■ complexity and anisotropy in the host population. The estimates for invasion threshold are linked to 
morphological characteristics of the hosts that can be used for determining the threshold for invasion in 

■ practical applications. 

o ■ 

1 Introduction 

•i-H , 

. One of the main questions in epidemiology regards the outbreak of epidemics, i.e. whether an infectious 

' disease can spread throughout a given ensemble of hosts or not. Many systems display a threshold for 

epidemics which divides the para meter space into regions where an outb reak may occur from regions 



epidemics wnicn aiviaes tne para meter space into regions wnere a,n outp reaK may occur irom regions 
where the disease cannot spread ( Murravl 2002 : Marro &: Dickman . 1999t ). Identification of the factors 



that determine such a threshold is of great importance in devising robust strategies to control the spread 
of disease through a population of susceptible hosts. Successful control deflects the system into the non- 
invasive region of parameter space so that a pathogen fails to invade. This simple concept can be applied 
to the spread of infection and disease, at a range of scales from the cellular, in which the host comprises a 
single cell, through populations of plants and animals in which a host equates with an individual organism, 
up to larger-scale systems, in which the unit of interest may be an individual field or farm for crop and 
livestock disease, a school, a village or other natural clustering for human disease. The challenge at each 
scale lies in dealing in a quantitative manner with factors such as stochasticity and heterogeneity, so that 
thresholds used to identify strategies for control are robust to these uncertainties. Stochasticity in disease 
spread is linked to the fact that susceptible hosts become infected only with a certain probability when 
challenged by inoculum from infected hosts under otherwise identical circumstances. Heterogeneity, on the 
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other hand, is associated with a range of factors that may differ amongst hosts. These include disorder 
in characteristics such as pathogen infectivity and host susccptibihty, and the spatial arrangement of 
hosts. For certain types of hosts, heterogeneity also includes the inherent morphological complexity (i.e. 
irregularity in shape) of the host. Morphological complexity is especially important in epidemiological 
spread where it affect the contact rate between contiguous hosts. Examples of morphologically complex 
hosts include dend r itic cells and plants, crop fields and farms, and individual clusters in social networks 



( Davis et ail . 120081 : iGonzalez eTai] . 120081 : ISoriano et oil . l2008t iBoender etal\ . l2007t lEisinger fc Thulkd . 



20081 ). 



While the problem of disease invasion has been extensively studied, both experimentally and theoreti- 
cally, most attention has been focused first on determinis t ic sys tems and increasingly on stochastic models 
( Murray . 2002 : Truscott fc Gilligan . 2003 : Gibson et all |200^ . Hete rogeneity has trad i tional l y received 



Levin fc DurrettI (Il996h; iNewmanI 



less a tt ention although there are some remarkable e x ceptions (e.g. . 

(l2002l) : ISander et all (l2002l |2003|) : [Cook et all (l2007l) : iBoender e/oll (|2007l ): lEisinger fc Thulkd (|2008l ): 
Milleil (|2007^ : lKenah fc RobinsI ()2007h : IVolzl (|2008l ) ). However, none of the previous work has established 



a link between heterogeneity and morphological features of the system leading to such a heterogeneity. 
In particular, the heterogeneity associated with the host morphology and the effects it may have on the 
features of the epidemics remain to be understood. 

In systems where the pathogen is transmitted between hosts due to their proximity, one can identify 
three main factors that determine the invasion threshold (see Fig.[T]): (i) spatial arrangement of the hosts 
in the population; (ii) morphology, i.e. the shape of the hosts, and (iii) the infection efficiency resulting 
from the net effect of the interplay between the pathogen infectivity and the host susceptibility upon 
contact. 

Sever al models have been pro posed for the dynamics of epidemics spreading by contacts between hosts 
(|Murravl . |i002; Liggett, 1985.: Marro fc Dickmanl . Il999t iHinrichsenl . I2OOO; Odor, 2004,). In such models, 
the hosts can be in different states, e.g. susceptible (S), infected (I), and recovered (or removed) (R) in 
the prototype SIR epidemiological model. The state of each host may change according to certain model- 
dependent rules. For instance, the SIR model assumes that infected hosts can infect others that are 
susceptible and then become recovered and fully immune to further pathogen attacks. In such models, 
the stochasticity and heterogeneity for the spreading process are treated in a simplistic way by using 
phenomenological probability densities for relevant parameters and thus not linking the host morphology 
and invasion threshold. 

In this paper, we establish a quantitative link between host morphology and the invasion threshold 
in an ensemble of hosts with realistically complex morphology. By analysing the conditions for epidemic 
outbreak in several systems with different degrees of configurational heterogeneity, we conclude that the 
invasion threshold is mainly determined by: (i) the average overlap between neighbouring hosts, (ii) the 
morphological complexity of hosts, and (iii) the host shape anisotropy. In particular, we demonstrate 
analytically and numerically that the resilience of the system to invasion increases with morphological 
complexity and anisotropy of hosts. This result is valid under very general conditions and is therefore 
applicable to a wide range of host ensembles. We show that irrespective of the degree of the host anisotropy 
the spreading process can be described in terms of a mean-field system, so that analytical estimates for 
the invasion threshold can be obtained. In addition, we complete our analysis by identifying several 
morphological characteristics that can be used for determining the threshold for invasion resulting in an 
epidemic outbreak. 



2 Methods 

We consider a set of N morphologically different branching structures, n = 1,2, .. .N, placed on all nodes 
i = 1, 2, . . . , of an L x L regular lattice with spacing a and nearest neighbour links only (see Table 
[??] for a summary of the notation used in the text). In particular, we deal with a triangular lattice 
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Figure 1. Example of a system formed by complex hosts represented, for concreteness, by planar 
neurons corresponding to the salamander retinal ganglion cells placed on the nodes of a triangular 
lattice with lattice spacing a in such a way that the somata coincide with the lattice nodes. The 
pathogen infests the surroundings (shaded area) of infected (I) hosts and, eventually, reaches the 
neighbouring susceptible hosts (e.g. amber susceptible (S) neuron on the right from I). The probability 
of infection of a susceptible host and the probability of a global epidemic outbreak depend on overlaps, 
J, between the infested region and susceptible hosts, and infection efficiency, k. The overlaps are 
dictated by the host morphology and spatial arrangement of hosts. The infection efficiency determines 
the effectiveness of the contact in terms of transmission of infection. Besides these factors, the inset 
shows that J, and thus the invasion threshold depend, in general, on the local orientation of the hosts, (p. 



(see Fig. [T]) because this arrangement corresponds to the organisation of hosts with the highest risk of 
epidemic spread, in the sense that any other 2D regular lattices with the same lattice spacing is more 
resistant to epidemic invasion (jlsichenkd . Il992l : IStauffer fc Aharonv , An additional advantage of 



the triangular lattice is that the next-nearest-neighbour links (ignored in our approach) are less likely in 
such a lattice as compared with other 2D lattices. A certain spatial configuration of generally anisotropic 
branching structures is fully defined in terms of the set of hosts placed on the lattice nodes, {ni\ (where 
rii is the host number placed on node z), and their orientations, {4>i]. The system of hosts can be either 
homogeneous, if morphologically identical hosts, rii = n, of the same orientation, 0^ = ^, are placed 
on the lattice, or heterogeneous if the hosts Ui and/or their orientations </),; are chosen, e.g. at random. 
Specifically, we consider three types of arrangements with different degrees of heterogeneity which allow 
the various factors contributing to the invasion threshold to be analysed separately: 

* Arrangements 1: both rii and (pi are distributed according to uniform distributions so that all 
possible hosts and orientations arc equally probable. These are highly disordered configurations. 

* Arrangements 2: the same host is placed on all the nodes with the same orientation (i.e., Ui ~ n 
and 4>i = (j) for all the nodes i). Such ordered arrangements highlight the role of the host anisotropy 
leading to different overlaps along different lattice directions. 

* Arrangements 3: the same host is placed on all the nodes and its orientation is drawn from a uniform 
distribution of width A(j) and mean value (f>. Such arrangements allow both the morphological 
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complexity and anisotropy of hosts to be analysed in a comparative manner and they contain 
arrangements of type 2 as particular cases with <j) = 4> and = 0. 

Table 1. Symbols and definitions used in the text. Indexes i and j span the L x L nodes in the 
(triangular) lattice. Index a spans the three main lattice directions, i.e., a = 1,2,3. Hosts are labeled 
by an index n ^ 1,2,..., N. The symbol pc stands for the bond percolation threshold. 



Host characteristics, arrangement, 


and overlaps 


p(r; 7i) 


Host density 
Host orientation 


a 

Jij ^ J p(r;ni)p(r - a;nj)dr 


Lattice spacing 

Overlap between hosts at nodes i and j 




Set of overlaps between all the pair of hosts 




Dispersion of J" associated with the shape complexity 


Dispersion of associated with the shape anisotropy 


Disease transmission 


k 


Infection efficiency 


T = 1 


Recovery time 




Transmission rate 


T,j = 1 - e-'3'^-^ = 1 - e-^-'^^ 


Transmissibility 


Invasion of infection 


Pinvik^ C^) 


Probability of invasion 


ac{k) = inf{a : Pinv(fc, a) = 0} 


Invasion threshold in terms of the lattice spacing 


fcc(a) = sup{fc : Pinv(fc,a) = 0} 


Invasion threshold in terms of the infection efficiency 


Invasion threshold in a mean-field system with overlaps 


Aki cx Vi 


(J) 

Non-mean-field contribution to k^ originated by the 
shape complexity 


Ak2 oc V2 


Non-mean-field contribution to k^ originated by the 
shape anisotropy 


Morphological characteristics relevant to the invasion threshold 




Radius of effective circular hosts in the mean-field ap- 


{^1,(^2,^3}, {di <d2< ds) 


proximation 

Lattice-adapted diameters 


SLAD = d2 


Second lattice-adapted diameter 



In order to investigate the spread of epidemics in such systems, we apply the dynamical rules of the 
SIR formalism. Infection and hence disease can be transmitted between infected and susceptible hosts 
with transmission rate /3, and infected hosts recover after a fixed time t. We assume that the value of t 
is time-independent, identical for all the hosts, and thus can be chosen as the time scale of the problem , 
T = 1 . Such homogeneity i n r provides a minimally safe bound for the invasion threshold (jKuulasmaal . 
1982 : Cox &: Durrett . 19881 ) (see brief explanation in Appendix E]). 

The transmission of infection from an infected host, rii at node i, to a susceptible nearest neighbour, 
rij at node j, separated by a unit-cell vector a, is a Poisson process with a transmission rate The 
value of Pij is assumed to be proportional to the overlap Jij between hosts i and j, i.e. /3ij = kJij, 
where k is the infection efficiency which accounts for the effectiveness of the overlap for transmission of 
infection. We shall set k to be identical for all the pairs of nearest neighbours. Possible variability in k 
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can be easily incorporated into the model but it does not change the main results qualitatively. 
The overlap between hosts is defined as 



Jy = J p{r; ni)p{r - a; ?ij)dr , (1) 



in terms of the host density, 



p(r;n,) = V <5(r-p) , (2) 



pSri 



where the position vector p scans all the points in the host structure n^. 

One of the main quantities involved in the SIR process is the transmissibility Tij defined as the 
probability that the pathogen is transmitted from an infected host at node i to infect a susceptible host 
located at node j during the life-tim e, r, of the infected node. For a Poisson process, the transmissibility 



is given by the following expression (jGrassbergeii Il983f ) , 



r„- = 1 - e-/^'^^ = 1 - e^'^'^-^ . (3) 

Therefore, the value of depends on the lattice spacing, a, which determines the overlap, Jy , and the 
infection efficiency, k. 

In a finite-size population, the SIR process lasts for finite time and after its termination two types 
of hosts, R and S, can be found in the system. The region of the recovered hosts can be mapped 
onto one of the clu sters in the bond p ercolation problem by mapping the transmissibilities {Tij} to the 
bond probabilities ( Grassbergeilll983l) . According to this mapping, the probability of invasion of disease 



starting from a single host, Pinv{k,a), is identified with the probability that this first infected host 
belongs to the infinite cluster of connected sites in the bond-percolation problem. The invasion threshold 
is formally defined as the boundary between the invasive and the non-invasive regimes characterised by 
Pinv > and Pinv = 0, respectively. This condition introduces a separatrix in the parameter space which 
can be given in terms of several (control) parameters of the system, such as, for example, the infection 
efficiency and lattice spacing. In terms of the infection efficiency, the invasion threshold is the value 
fee = sup{fc : Pinv = 0} such that relatively small values of fc < fcc correspond to non-invasive regime (i.e., 
-Pinv = 0); while larger values, k > k^, describe the invasive domain (i.e., Pjnv > 0). Similarly, the critical 
lattice spacing, Oc = inf{a : Pinv = 0}, splits the range of lattice spacings into two regions: a < ac and 
a > Oc corresponding to invasive and non-invasive regimes, respectively. The host morphology affects 
implicitly (through the overlaps between hosts {Jij}) both fcc(a) and ac{k). The morphological variability 
of the hosts together with the heterogeneity in their arrangement makes the set { J^ } disperse in general. 
The overlap can then be regarded as a random variable J taking the values {Jij}- We will describe J 
in terms of its average, {J), and deviations from (j7} originated by its dispersion. Consequently, it is 
convenient to split the expression for invasion threshold into two contributions: 

fc, = kl + Afc, (4) 

where and Afc are associated with the average and the dispersion of J ^ respectively. For the spatial 
arrangements of branching hosts considered below, the source of dispersion in the overlaps is two-fold: (i) 
Vi , which is due to variability arising from complexity in the shapes of different hosts for arrangements 
of type 1 and (ii) V2, which is due to variability in shape anisotropy for arrangements of type 2. Corre- 
spondingly, the value of Afc can be split into two components, i.e. Afc = Afci -I- Afc2, where Afci oc Vi and 
Afc2 oc V2. Both complexity and anisotropy contribute to dispersion of the overlaps for arrangements of 
type 3. 

For numerical illustration and concreteness, in this paper we use a set of N (N = 51) neurons (Fig.[T|) 
corresponding to the salamander retinal ganglion cells ( Ascoli . 20061 ). which are mostly planar, as typical 



representatives of complex branching structures. In this case, the vector p introduced in Eq. ^ scans all 
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the pixels defining the digital image of each neuron. Technically, the (5-functions in Eq. ^ are replaced 
by Gaussians of width comparable to the pixel size. This broadening mimics the (diffusive) spreading of 
the pathogen around the host. While we are not aware of documented examples of pathogen transmission 
in these structu res, they are represe n tative of a broad class of s tr uctures that are known to transm i t virus 



infections (e.g. iLaVail et al\ (|1997^ : lEhrengruber et al\ (|2002D : IChen et al\ (|2007f ): ISamuel et al\ (|2007^ 



). The importance here is to use the published data on complex morphology to test the general methods 
introduced below. 



3 Results 

In this section, we analyse the invasion threshold resulting in each of the spatial arrangements listed 
above and propose morphological characteristics for description of the invasion threshold. 



3.1 Arrangements 1. rij and random 

The invasion probability in host arrangements of type 1 is presented in Fig. [2] In particular, the inset 
shows the dependence of the probability of invasion on the infection efficiency. The value of Pinv(fc) is zero 
at small values of k and becomes positive above the threshold, k > kc{a). As expected, fcc(«) increases 
with the distance between hosts (cf. the curves marked by different symbols in the inset), a, since the 
infection mechanism between neighbouring hosts should be more efficient in order to invade the system 
with larger lattice spacing. However, the curves for Pi^v at different a collapse onto a single master curve 
(see Fig. [5]), if plotted as a function of an alternative control parameter, the aver age transmiss i bility . 



(T), which plays the same role as bond probability in the pe rcolation prob l em ([Sander et all 12002 , 



2003 ). This statement was first demonstrated heuristically bv ISander et al. ( 20021 ) for SIR processes 
with heterogeneous transmissibility. Indeed, it is easy to go beyond heuristic arguments and prove that 
the statement is valid if the tran smissibilities I Till are independent and self- averaging qu antities (see 
demonstration in Appendix E and |Kuulasm3 (|l982[ ): iMille^ (|2007t ): IXenah fc Robinl (|2007t ) for cases in 
which the transmissibilities are not independent and such description is not appropriate). The collapse 
in Fig. [2] suggests that transmissibilities in arrangements of type 1 satisfy these conditions. In fact, the 
master curve defined by the collapse coincides with that obtained for the effective homogeneous 'mean- 
field' system with transmissibility (T) between all the nearest neighbours (see the dashed curve in Fig. [2] 
which is identical for all lattice spacings). This means that the heterogeneous system is equivalent to a 
homogeneous 'mean-field' one for which the invasion threshold can be determined in an efficient manner 
by solving the following equation, 

a)=Pc, (5) 

0.347 is the bond perco l ation threshold in an infinite (i.e., L oo) triangular lattice 



where 



2lL 



(jlsichenkd . 119921 : IStauffer fc Aharonvl Il994f). The value of Pr. provides the minimally safe bound for 



invasion threshold in other 2D lattices ( Isichenkol 1992 : Stauffer &: Aharonvl . Il99-i ). Since the average 
(T) converges very fast with L to its limiting value for L ^ cx), the value of kc{a) estimated from Eq. ([5]) is 
representative for macroscopic systems. The dependence of infection efficiency on lattice spacing defines 
the phase boundary in the (fc, a) plane (see solid circles in Fig.[3l^a)) between the invasive and non-invasive 
regimes. The threshold ac(fc) provides the same separatrix. 

The influence of the host morphology on the invasion threshold can be better understood by analysing 
the dependence of kc on the overlaps J . The first term in Eq. ([4]) can be found by solving Eq. ([5]) for 
the homogeneous system in which all the overlaps are replaced by its mean value. {J), i.e. 



|ln(l-Pc)| 
{J) 



(6) 
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Figure 2. Probability of invasion, Pinv, for morphologically complex hosts in disordered arrangements 
of type 1 on a lattice of size L x L = 200 x 200. The main plot displays Pinv as a function of the average 
transmissibility, (T), for heterogeneous systems with different lattice spacings a (marked by different 
symbols) and for a mean- field system with homogeneous transmissibility (T) (dashed line). The 
bond-percolation critical probability, pc, marked by arrow gives the invasion threshold in the 
thermodynamic limit. The inset shows the invasion probability as a function of the infection efficiency k 
for different values of lattice spacing marked by the same symbols as in the main figure. 



It can be rigorously shown that the value of fcj? underestimates the threshold, i.e. kc > fcj? (cf. circles and 
dashed line in Fig. [31(b)). To prove this inequality we write the expression ([5]) in terms of the overlap as 
1 — (e"'^"'^) = Pc- Applying then the general inequality {e~^'"-^) > e~'^=^'^^ to the above relation we obtain 
kc > — ln(l —pc)/{J) which reduces to kc > k^ after using the expression ([6]) for k^. This inequality 
implies that the contribution Afc to kc associated with the dispersion of the overlaps makes systems more 
resilient to epidemic invasion. An approximate solution of Eq. ([Sj obtained by keeping the first correction 
to fc° gives the low-bound estimate Afci of A/c, i.e. 



which is proportional to the variance of the overlaps, Vi = {J"^) — {J)"^ ^ with other central moments of 
higher order being dropped (see the derivation of Eq. ([7|) in Appendix ID|) . The upper inset in Fig. [3] 
shows that the value of Vi is non-zero, as expected from the comparison of kc and fcj? plotted in the main 
figure. 

The low-bound estimate of Afc given by Eq. ([7]) provides a safe threshold, fcj? 4- Afci, for the actual 
system. We expect this to be the case for a wide class of morphologically complex hosts so that describing 
J in terms of its average and standard deviation provides a minimally safe bound to the invasion threshold 
in heterogeneous arrangements. 

3.2 Arrangements 2. Identical and 0j at all the nodes 

In ordered arrangements of type 2, all the hosts are identical branching structures with the same ori- 
entation. The host overlaps along the main lattice directions, {Ji, J2, J3} (see Figs. IHa)), depend on 
the lattice spacing a, the host used as a motif, n, and its orientation, (p. For a given n, the set of all 
possible overlaps obtained by varying a and (j) within their respective domains define the host overlap 



Afc > Afci 



{k^cY 



(7) 
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Figure 3. Invasion tiiresliold for tiie system of brancliing liosts. (a) Representation in terms of the 
infection efficiency, fc, and the lattice spacing, a. The threshold kc{a) for arrangements of type 1 is 
shown by circles. The shaded region corresponds to the statistically possible values for the estimation of 
kc{a) in terms of effective circular hosts with homogeneous overlaps defined in Eq. ([H)) . Diamonds 
indicate the average thresholds kc{a) (solid symbols) and adk) (open symbols) for arrangements of type 
3. (b) Invasion threshold fee as a function of the average overlap (jj) for arrangements of type 1 
(circles), type 2 with orientation = (diamonds) and type 3 with mean orientation = and two 
widths of uniform distribution, A0 = 1 (solid squares) and At/) = 27r (open squares). The branching 
host used as a motif for arrangements of type 2 and 3 is displayed in the figure. The dashed line 
represents the dependence of vs (J) given by Eq. (|6]). The insets show the dispersions Vi and V2 of 
the overlaps associated with the disorder and anisotropy, respectively, corresponding to the same 
arrangements as in the main figure (the symbol code is the same as in the main figure). 
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locus. Fig.Ul^b) shows the overlap locus corresponding to a typical branching host displayed in Fig. |31Ja). 
Each configuration with given a and 4> is represented by a point (Ji(a, (/)), J2(a, 0), J3(a, (/>)) belonging 
to the overlap locus. The inherent anisotropy of the host is reflected in the dispersion of the overlaps, 
Ji 7^ J2 7^ J3 7^ Ji, and gives rise to a significant deviation of the overlap locus (surface in blue) from 
the straight line (in black) corresponding to the locus for isotropic host ( Ji = J2 = J3). 

Similarly to arrangements of type 1 , the spread of disease in ordered arrangements can be mapped onto 
the bond percolation problem but now with anisotropic bond probabilities corresponding to the values of 
transmissibilities, Ti 7^ T2 7^ T3 7^ Ti, along the main lattice directions defined as Ti = 1 — e"*'"^' for i 



1, 2, 3 . From this mapping, the invasion threshold is defined by the following condition (|Svkes &: Essam 



19641) 



I-T1-T2-T3 + T1T2T3 = . (8) 

This condition can be recast in terms of the infection efficiency and overlaps, as g( Ji, J2, J3, fcc) — 0, 
where 

5(Ji, J2, J3, fc) = 1 - e-"('^i+^^) - e~''^''+-'-'^ - e-*^('^^+^^) 

^^g-fc(./l + ./2 + ./3)_ 

For given value of k, the condition .g = defines a critical surface (in red in Fig. IH^b)) in the overlap 
space. Therefore, the intersection of the critical surface with the overlap locus (in blue) defines the 
invasion threshold, ac(fc, 0). On the other hand, for a given value of a and 0, the critical infection 
efficiency kc{a^4}) is given by the value of k which generates a critical surface containing the point 
(Ji(a, 0), J2(a, (^), J3(a, 0)). 

In practice, it is useful to consider the minimally resilient thresholds fc™'"(a) ~ mm^{kc{a, (p)} or 
a™^^{k) = maxcf,{ac{k, (f>)} which ensure that the system is safe for any orientation if k < fc"""(a) or 
a > a^^^{k), respectively. The invasion thresholds are host-dependent and define different separatrices 
in the (a, k) plane for each host. Fig. EJa) shows the average thresholds A:™'"(a) and a™^^(fc) over all the 
branching hosts {n}. The phase boundary, a™^^(fc), gives the safest estimate for the invasion threshold, 
which is a consequence of the multivalued nature of adk, (p) in contrast to single- valued function kc{a^ (f) 
(see Appendix [Cl for more detail). 

The effects of the host anisotropy on the invasion threshold can be analysed in a similar way as for 
configurations of type 1 by investigating the dependence of kc on the overlaps. Similarly, the critical 
infection efficiency, is given by Eq. (|4]) with the mean-field value fcj? evaluated for the system with mean 
overlap {J) = (Ji + J2 + Jz)/^- The dispersion in the overlaps results in the following approximate 
expression for Afc (see detailed derivation in Appendix [D|) . 

Afc > Afc2 = \^V2, (10) 

where the quantity V2 = z(i+p ) Eq=i {^a — {J))^ accounts for the anisotropy of the overlaps. The value 
of Afc2 is non- negative meaning that kc > fcj?, i.e. the anisotropy in host shape makes the system more 
resilient as compared with the system of isotropic hosts with the same mean overlap. 

As an example. Fig. El^b) shows that kc > k^ in the system of branching hosts, in agreement with the 
predicted behaviour (cf. diamonds and the dashed line). The inset shows the corresponding dispersion 
V2. Given that the deviation of kc from fcj? is small in the system of branching hosts considered, the 
correction Afc2 is in fact a good approximation to the actual deviation Afc. This is the expected behaviour 
for systems of hosts with moderate anisotropy. 

3.3 Arrangements 3. Identical and random 0j 

In the two previous sections, it has been shown that both the disorder and anisotropy of the hosts 
make systems of branching hosts more resilient against epidemics. The arrangements of type 1 and 2 are 
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Figure 4. Arrangements of type 2 (an identical host witli the same orientation placed at all the 
nodes), (a) Unit ceU of two different configurations constructed by using the same branching host. The 
orientation, 0, and spacing, a, corresponding to each case are indicated in the schematic reference 
frames displayed at the top. As shown in the left frame, there are three different values of the overlaps, 
{ Ji, J2, J3}, corresponding to each of the main directions in the lattice, (b) Space of overlaps, 
{Ji, J2, Js)- Each configuration with different a and (j) is mapped into a point in this space. The set of 
overlaps corresponding to all the possible configurations obtained for a given host defines its overlap 
locus (surface in blue). The deviation of the blue surface from the straight black line representing 
overlaps between isotropic hosts { Ji — J2 = J3) shows the degree of anisotropy in the overlaps for a 
typical branching host (shown in (a)). High values of the overlaps correspond to small lattice spacings 
a. In particular, the points labeled as CI and C2 correspond to the configurations shown in (a). The 
critical surface defined by g(Ji, J2, J3, k) = (see Eq. ^) is shown in red for two values of the infection 
efficiency: k — 0.02 and k = 0.008. For a given value of k, the points in the overlap locus below/above 
the critical surface correspond to safe/ vulnerable configurations. The intersection of the overlap locus 
with the critical surfaces parameterized by k determines the critical threshold ac{k, 4>). The critical 
infection efficiency fee (a, (f) corresponding to a configuration with lattice spacing a and orientation is 
given by the value of k related to the critical surface containing the configuration point 
{ Ji(a, (/)), J2(a, </>), J3(a, </))}. 
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extreme cases in the sense that the first type highlights the morphological complexity {Vi > and V2 = 0) 
whereas the second type highlights the effects of the anisotropy of the hosts (Vi = and V2 > 0). In 
generic arrangements, such as those of the type 3 defined above, the two effects are present. By assuming 
that the transmissibilities between different pairs of hosts are independent, it is possible to show that the 
behaviour of the actual heterogeneous system is equivalent to that of a mean-field homogeneous system 
with anisotropic transmissibilities {(7i), (72), (Ts)}, where {Ta) is the average of the transmissibility along 
the direction a in the lattice. A proof for this statement, which represents a generalisati o n to system s 



3 syi 

with anisotropic transmissibilities of the mean- field description suggested in lSander et al\ (|2002i |2003[ ). 
is given in Appendix |B] The equation for the invasion threshold in this case is 

1 - (Ti) - {T2) - {%) + {Ti){T2){%) = 0, (11) 

which generalises the formulae ^ and ([8]) in such a way that Eq. ^ corresponds to the particular case 
of Eq. (fTTj) when (T) is the same along all the directions and Eq. ([8]) emerges when there is no disorder 
in the anisotropic transmissibilities. 

As for the previous arrangements, the critical infection efficiency obeys Eq. ^ with both heterogeneity 
and anisotropy contributing to Afc. The lowest-order approximation to the lower bound for Afc is given 
by the relation Afc > Afcs = Afci + Afc2 (see details in Appendix |D|) where Afci and Afc2 are defined in 
Eqs. ([7]) and (fTO|). respectively, with the dispersions terms generalised to 

yi=\jl{{Ji)-{Jc.?), (12) 
^^-3(TT^t(a)-(^))^- (13) 

^ ' a— 1 

Fig. [3l^b) shows that the increase in orientational variability A(/) brings additional heterogeneity in the 
system and thus results in a decrease of V2 but this does not necessarily induce a decrease in the value of 
critical infection efficiency since the contribution Vi may increase. This illustrates the interplay between 
the role of disorder and host anisotropy. 



3.4 Description of the invasion threshold in terms of morphological charac- 
teristics 

In the previous sections, we have established a link between the invasion threshold and host overlaps 
characterised by the first moment, {J^), and the deviations from the mean, Vi and V2. The missing link 
between the host morphology and invasion threshold can be recovered by studying how the morphology 
affects overlaps. Here we show that both {J) and the anisotropy of J (i.e. V2) can be well described 
in terms of a reduced number of morphological characteristics. In contrast, a proper description of the 
disorder-induced contribution from Afci in terms of a reasonably small set of morphological characteristics 
is hardly possible and requires instead knowledge of the spatial host density. p(r) (sec Eq. However, 
we can ignore the disorder-induced contributions from Afci in order to obtain a safe lower bound for 
critical infection efficiency and thus connect this quantity with morphological characteristics of hosts. 

We start the analysis by considering arrangements of type 1. The overlaps arc statistically isotropic 
(i.e., V2 — 0, Fig. |3l[b)) so that the invasion threshold can be described in terms of a mean-field system 
with homogeneous and isotropic overlaps, and corresponding transmissibilities (T). The mean- field 
system consists of effective circular hosts of radius r™^ with overlap J'™^. These effective circles are fully 
described by the radius, r™^, and density, p™'(r), which is positive for < r < r™^ and is zero otherwise. 
The functional form of p™^(r) is the same as of real branching hosts in [0,r™^], i.e. p{r) oc r'^i~'^ where 
df stands for the fractal dimension (see Appendix [El for more detail). The radius of the effective circles 
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is defined as r™^ = ^^jf-^ (r) to ensure that the average radius of the effective circles coincides with the 
average radius of the actual hosts, (r) (see Appendix |E|) . Under these assumptions, the expression for 
the overlap between neighbouring effective circles is given by 



valid for small overlaps when 2r'°^^/a — l <C 1. The parameter AI is a normalisation constant (see Appendix 
IE)1 . The substitution of this expression for {J') into Eq. ^ gives the required link between the invasion 
threshold and the effective radius, k™^ oc (r™^/a)^/^(2r™Ya — In fact, for the particular system of 

branching hosts studied here, this is a very good estimate for the value of the critical infection efficiency, 
kc ^ k™^, as demonstrated in Fig. [3] (the solid circled line falls inside the shaded gray area representing 
the set of possible values for fc™^ within statistical errors). This means that for the arrangements of type 
1, the mean-field homogeneous system of effective circles can be reliably used for estimating the phase 
boundaries. 

In the case of ordered (type 2) and partially ordered (type 3) arrangements, the overlaps along distinct 
lattice directions can be significantly different and thus crucial for evaluation of the invasion threshold 
in contrast to the mean overlap which is important for disordered (mean-field like) arrangements of type 
1. Therefore, instead of a single characteristic such as mean radius of the effective circles for disordered 
arrangements, we introduce an ordered set of linear sizes of branching hosts along the lattice directions, 
i.e the set of lattice-adapted diameters, {^1,^2,^3} (di < d-i < d^; see illustration in Fig. [5]). The second 
lattice-adapted diameter (SLAD), d2, plays the most important role for finding the critical value of, e.g. 
lattice spacing, and thus estimating the invasion threshold. This is due to the fact that Eqs. ([5]) and 
(ITT|) have a solution for the critical threshold only if the overlaps at least along two directions are finite. 
The overlaps become greater than zero if the lattice-adapted diameters are comparable or greater than 
the lattice spacing and thus solution of Eqs. ^ and PT|) exists if a < d2. This qualitative analysis 
suggests the existence of strong correlations between the SLAD and Uc- Indeed, we have observed such 
correlations, i.e. Oc — d2, for ordered arrangements of type 2 (see Fig. [5|). The description of Oc in 
terms of SLAD gets worse for small values of k (see the black squares in Fig. O when the strong overlaps 
between hosts should be achieved at criticality and thus the interior density of the hosts, rather than the 
diameter only, becomes important. 

4 Discussion and Conclusions 

Using the framework of an SIR epidemiological model, we have investigated the transmission of infection 
and spread of disease in systems of hosts with realistically complex morphology. Our main finding is 
that the greater the irregularity in the host morphology, the more resilient is the population to epidemic 
invasion under otherwise identical conditions (for instance, identical spatial arrangement of hosts). We 
derive a safe lower bound for the invasion threshold, which has been obtained for branching hosts with 
independent transmission rates placed on a triangular lattice. We have shown mathematically that this 
bound holds for all other 2-D topological arrangements of hosts with nearest-neighbour transmission and 
even in the case when the transmission rates between different hosts are correlated with each other. In 
particular, irregularity in the host positions, i.e. small random displacements of hosts from lattice nodes 
bring correlations in transmissibilities and thus make the system more resilient. Of course, for some 
real systems, the assumption about nearest-neighbor transmission may be violated by the presence of 
short-cuts between remote nodes due, for example, to wind or animal motion. I n this case, it is known 



that the system becomes less resilient to epidemic invasion (jSander et all 120021 ) and the bounds given 
above are no longer valid. 

We have used a set of planar neurons to illustrate the effects of complex branching hosts on the spread 
of infection. We have identified two sources of heterogeneity in the systems considered: (i) morphological 
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Figure 5. Dependence of the critical lattice spacing, a™^^{k) on the SLAD, d2, for all the branching 
hosts in arrangements of type 2. Different symbols refer to different values of k as marked and each 
point for a particular symbol corresponds to an individual branching host. The solid lines represent the 
linear regression fit for each value of k (e.g. ac{k) = 14.1 + 0.9Ad2 with correlation coefficient ~ 0.97 for 
k = 1). The inset defines graphically the lattice-adapted diameters {c?i, (i2, rfs}- 



complexity of hosts and (ii) the host shape anisotropy. Both contribute to the resilience of the system 
against epidemic invasion and can be described by means of two morphological characteristics of hosts, i.e. 
by the mean effective radius and the second lattice-adapted diameter. Such characterisation is not exact 
in general but it provides a safe bound to the invasion threshold. The main conclusions are generic and 
remain valid for any type of morphologically complex hosts. In particular, the methodology introduced 
here and the bounds for resilience to invasion apply to the transmission of infection in other 2-D systems 
with analogous disorder expressed in the host morphology and anisotropy. Examples include the spread of 
plant disease through contacts between adjacent plants in a field, orchard or forest in which the host plants 
frequently occur on a 2-D lattice. Here the contact structures between nearest-n eighbours are d etermined 
by overlap of shoots (for aerial pathogens) or roots (for soil-borne pathogens) (jGilliganl . |2008| ) in which 
the 3-D structure of the plant can be collapsed onto a 2-D framework when considering transmission of 
infection between nearest-neighbours. In principle, a similar analysis could be performed for ensembles 
of morphologically complex 3-D hosts arranged on a 2-D lattice that takes explicit account of the three 
dimensional host structure. Expressions for the invasion threshold are not known analytically in this case 
and the precise definition of quantities directly linked to host morphologies such as fractal dimensions 
or average radius are system-dependent that require further study. Other potential applications of the 
methods include analysis of the transmission of infection via the 'morphology' of contacts between clusters 
of susceptible hosts in social networks, in which the clusters can be approximated by a 2-D lattice (cf 
recent work on percolati on models for the spread of plague through gerbil populations in lattices of 
interconnecting burrows ( Davis et al. . 20081 )). 



The present work opens several possible directions of further research for practical applications in 
considering disease control strategies and for basic understanding of epidemic spread involving heteroge- 
neous transmission of infection. For instance, our analysis suggests new ways for control of epidemics in 
real systems where host morphology is inherently complex. For example, in a system where an epidemic 
is active, a treatment enhancing the anisotropy in the transmission rates, would be more efficient as 
compared with reduction of transmission rates in all directions, i.e. isotropically. Such might arise with 
the deployment of microbiological biological control agents to restrict the spread of infection of soil-borne 
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pathogens in plant pop ulations ( Gibson et al. , 19991) . Biological control agents often exhibit marked vari- 



ability in performance (jGibson et all Il999l ). faiUng to provide isotropic control of infection on targeted 
hosts, for example due to uneven colonization of roots by microbial antagonists deployed as biological 
control agents. Our results suggest that increasing the degree of anisotropy due to these organisms may 
yet contribute to success in controlling invasion. Anisotropy may also be fostered in social or animal 
systems by preferential treatment of some, rather than all, connections between clusters. The fin dings in 



the cu rrent manuscript are also important for better understanding virus tracing of neurons (e.g. iLoewv 



(|l998h '>. which is a biological staining method where virus propagation from neuron to neuron is used as 



a means to histologicaly mark the interconnections, so that they become visible to the microscope. More 
specifically, our results imply that virus tracing might be not so effective in the case of anysotropic or not 
so complex neuronal cells, which could therefore be overlooked by this type of marking. Analogue effects 
can be also important in transneuronal spreading of virus in order to deliver gene therapy (,Oztas, 2003). 

Finally, the approach presented here is relevant to epidemics for which an SIR model is suitable and 
a mapping to ordinary percolation exists. While many diseases can be described by the SIR framework, 
others cannot. It would be interesting to analyze the effect of host morphology within the framework 
of a different family of epidemiological models with final state not immune to the disease (such as the 
susceptibleinfectedsusceptible model). Such an extension is challenging since the mapping to ordinary 



susceptiblemiectedsusceptible model j. bucn an extension is cnallcngmg since tne mapping to ordinary 
percolation is no longer possible requiring tools such as directed percolation ( Marro fc Dickmanl . Il999t 
Hinrichsen . I2OOOI) from non-equilibrium physics. 



Appendices 

A Homogeneous recovery times 

In the present work, we have assumed that the recovery times of the hosts are homogeneous. This 
assumption provides a minimally safe bound for the invasion threshold for epidemic outbreak so that the 
systems with disorder in r are safer i.e. less likely to be invaded. Indeed, the probabil ity of invasion P--^ 
in an SIR process ^ with heterogeneous t and average transmissibility (T) satisfies ( Kuulasmaal 1982t 



Cox fc Durrettl . 119881 ): 

p.-^ < p.^ (^ 

^ inv — ^ invi \^'~'/ 

where Pj^^ is the probability of invasion in an SIR process B with the same average transmissibility as in A 
but with homogeneous r. The above inequality implies that the non-invasive region for the heterogeneous 
(in r) system is wider than for the homogeneous one. 

B Validity of the mean-field description of heterogeneous sys- 
tems 

The validity of the mean- field description leading to the t hreshold cond i tion ( 2) in the main text in systems 
with isotropic disorder was heuristically suggested by ISander et al. (l2002h . In fact, it can be proven 



more generally that the mean-field description is valid both in the presence of isotropic and anisotropic 
uncorrelated disorder provided the overlaps between different pairs of hosts (and thus transmissihilities) 
are self- averaging independent quantitie)^. In order to prove this, let us interpret Ty as the conditional 
probability, Tij{J) = Pii — > j\J), for infection to be transmitted from host i to host j given the overlap 
J which is a random variable. In the anisotropic case, the overlaps along different lattice directions are 
distributed according to distinct and independent probability densities, e.g., {/a(J),Q! = 1,2,3} for a 



physical characteristic of a disordered system is said to be self-avera ging if its average over several configurations of 
disorder coincides with its average in a single infinitely large configuration l|Sornett'3 , 12OOOI ) 
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triangular lattice. The probability that the disease is transmitted from an infeeted host at node i to one 
of its neighbors at node j along the direction a is given then by the following expression, 

(t}-^)= jT,,{J)UJ)dJ (16) 

which is identical for all the pairs along direction a and thus coincides with the average trans- 

missibility (7^) = {T.j;"^) if the transmissibility is a self-averaging quantity. We have checked that the 
self-averaging condition indeed holds for branching structures by demonstrating that the probability of 
invasion does not depend on the particular configuration of hosts. This is illustrated in Fig. |6] where we 
show that the invasion probability curves (Pinv vs (7~)) collapse for different realizations of disorder. 




As a consequence of the above analysis, the system can be mapped onto an anisotropic but homoge- 
neous mean-field system with transmissibilities, {(7^)}. Therefore, Eq. (8, main text) giving the invasion 
threshold for an anisotropic disordered system can be obtained from Eq. (5, main text) by replacing 
the anisotropic transmissibilities {Ta} by the mean values {(7^)}. In the particular case of isotropic 
triangular system, (Ta) ~ (T), so that Eq. (8, mai n text) reduces to 1 — 3(T) -|- (T)^ = whi ch has the 
solution (T) = Pc, where Pc ^ 2 sin(7r/18) ~ 0.347 (jlsichenkol Il992l: IStauffer fc Aharonvi [l993) . 



C Multivalued behaviour of critical lattice spacing 

The invasion threshold in the ordered arrangements of type 2 (the same host is placed on all the nodes 
with the same orientation, </)) is given by the condition g = (see Eq. (6) in the main text). The 
function g depends on the overlaps along the main lattice directions, {Jaia,4>),a = 1,2,3} . For fixed 
orientation </>, the overlap Ja(a, </>) between branching structures is generally a non-monotonic function of 
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the lattice spacing (see blue and red curves in the inset in Fig. [7])- In this figure, we show that the non- 
monotonic character of the overlaps results in non-monotonic behaviour of the critical infection efficiency, 
fcc(a, (/)), and confers a multivaluated dependence of the critical lattice spacing, ac(fc, 0) versus k. As a 
consequence, the minimally resilient threshold in terms of the lattice spacing, a™^'^(fc) = max0{ac(fc, (f))} 
(see the illustration for a^^^{k) in Fig. [7] where this quantity is represented by the red curve) gives 
a safer boundary as compared to the minimally resilient threshold in terms of the infection efficiency, 
A:™'"(a) ~ mmcf,{kc{a, (j))} (the black curve), i.e. the red curve in Fig. [7] either coincides with or is to the 
right of the black one. 




Figure 7. Infection efficiency versus lattice spacing for fixed orientation of branching hosts represented 
by a planar neuron shown in the inset. The definition of the maximally resilient threshold a^^^(k) (red 
curve) is illustrated for a particular value of k for which ac(fc, (/>) takes three different values. The inset 
shows the dependence Ji(a) along one of the lattice directions for the arrangement of type 2 obtained by 
using the displayed host as a motif. The shaded (brown) region indicates the overlaps for all the possible 
orientations, 0, of the host. Blue and red lines correspond to Ji(a) for ^ = and </> = 1, respectively. 



D Analytical estimates for the invasion threshold 

In this Section, we derive an approximate analytical expression for kc for arrangements of type 3 and 
extend this result to other arrangements. 

We start from Eq. (8) in the main text which is valid for arrangements of type 3. Expressing the 
transmissibilities, {7^, a = 1,2,3}, in terms of anisotropic overlaps, {Ja, ol = 1,2,3}, gives 

1 - {e-^^^^)(^e-^^-^^) - {e-^-^^){e-^':^^) - {e-^-^^){e-''''^^) + (e-'=<='^0(e"'''=-^')(e~'''=-^^) = . (17) 

The overlaps, J7q = {J^a) + 5Ja, along lattice direction a are random values characterized by the mean 
value {Jo) and deviations from the mean, 5J^a- Anisotropy brings an additional source of dispersion, 
AJTq = {Jo) — {J), where the overall mean overlap is {J) = \ a^ii'^a) ■ Expanding Eq. (|17p in small 
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^Ja *C {J a) fii'st and then in AJTq <C and keeping the leading eorreetions to the the mean- field 
value, resuhs in the foUowing expression: 

fc2g-2fc.(J> ^0 _ (18) 

In order to derive an analytical estimate for the critical infection efficiency from Eq. (jlSp . it is con- 
venient to separate the mean- field contribution, fc^, i.e. fee = fcj? -|- Afc, and evaluate the value of Afc by 
expanding Eq. (jlSp in Afc ^ fcj?. The zero-order term gives the condition for the critical mean- field value 
of the infection efficiency, 

l_3e-2fc«U>+e-3/c«U) ^0 , (19) 

which has the solution fc-" = — '"^[^^^^^ , where = 2sin(7r/18) ~ 0.347. 
The next-order terms give an estimate for Afc, 

Afc>Afc3==||^(Fi+l^2) , (20) 

where V\ and Vi are defined in the main text by Eqs. (9) and (10). 

The estimates of Afci and Afc2 for arrangements of type 1 and 2, respectively, can be obtained in a 
similar way as particular cases of the derivation given above. 



— fe {J) \ 1 



E Mean-field description of complex shapes 

In this section, we show how the complex shapes can be approximately represented by effective regular 
shapes, i.e. effective circles. 

In disordered arrangements of type 1, an arbitrary host n,; from the set of nodes {n\ (n ~ 1, • . ■ , 
is placed at each lattice node i at uniformly random orientation (/>,;. The invasion threshold can then be 
described in terms of a mean-field system in which all the hosts have the same density, (p) (r) , where 

defined as the average of the host density p(r; n) over the iV morphologically different hosts in the 
ensemble and their orientation. Such averaging eliminates the dependence on the orientation so that 
(p)(r) depends on the radial distance, r — |r|, only. By definition, the integral of the density over the 2D 
plane gives the average "mass" of hosts, i.e. the mass of the effective host in the mean-field system, 

M = / (p)dV . (22) 

The mean radius of the effective host is then given by 

(r) = / r(p)dV . (23) 

The mean density can be easily obtained by numerical integration of Eq. (PT|) . Fig. [H] shows that the 
mean density for branching structures appears to decay with radius according to a power law, (p) ^ r^f^"^ 
with df = 1.5 ± 0.1 (truncated at large values of r by an exponential cut-off which accounts for the 
finite radius of the neurons used here for illustration). This is the expected behaviour for branching 
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Figure 8. Average density for brandling structures for disordered arrangements of type 1 versus r. The 
dashed line shows a pure power-law behaviour with exponent df — 2 = —0.5. The inset illustrates the 
overlap of two effective circles of radius, r™^, placed at a distance a apart from each other. The overlap 

'o];r e [r,„i„(^),r"^f]}. 



region is defined by the following conditions in polar coordinates (r, 6): {0 G [0q, —6[ 



structures which are non-compact objects and are then characterized by a fractal dimension df smaller 
than dimensionality d = 2 of the embedding space. 

The exponential decay of the mean density at large values of r can be ignored in the first approximation 
so that the hosts in the mean-field description can be represented by effective circles of finite radius r™^ 
with density, p™^(r), given by 

p-f(,)^P^ 2. ^ ' (24) 
[0, r > r"'' . 

The mean-field radius is defined by the following condition, 

27r r'"' 



(r) = d9 rp^Hr) rdr = '—r'^' . (25) 

J J df + 1 



The value of the mean radius (r) in Eq. (|25|) can be evaluated for the heterogeneous system using Eq. ([23|) 
and thus used for definition of r™^ to ensure that the mean radius of the effective circle coincides with 
the mean radius of actual hosts. 

The overlap between two effective hosts represented by mean-field circles placed at a distance a apart 
from each other (see inset in Fig. is given by Eq. (12) in the main text, i.c by the integral, 

'^"^'-(^^^^) / / r'^'-'{^r^+a--'2arco.ef'\r^ (26) 
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where smd r,„in arc defined in Fig. [S] If the overlap is small, the above integral can be approximately 
evaluated by performing an expansion in powers of the linear overlap S = 2r™^/ a — 1 ^ 1 , giving 



In the system of branching hosts used for illustration in the main text, M = 1282.4 stands for the average 
number of pixels in the digital images, df = 1.5 ± 0.1, and (r) = 50 ± 5. The effective radius is then 
r™f = 83 ± 10 and the approximation to the overlap becomes J'°f ~ (0.024 ± 0.008)(53/2. 
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